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Abstract The transition to turbulence in flows where the laminar profile is linearly stable requires 
perturbations of Unite amplitude. “Optimal” perturbations are distinguished as extrema of certain 
functionals, and different functionals give different optima. We here discuss the phase space structure 
of a 2-d simplified model of the transition to turbulence and discuss optimal perturbations with 
respect to three criteria: energy of the initial condition, energy dissipation of the initial condition 
and amplitude of noise in a stochastic transition. We And that the states that trigger the transition 
are different in the three cases, but show the same scaling with Reynolds number. 
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INTRODUCTION 


In parallel shear flows like pipe flow, plane Couette 
flow or Poiseuille flow and in boundary layers like the 
asymptotic suction boundary layer or the Blasius profile, 
turbulence appears when the laminar profile is linearly 
stable against perturbations [ij. Accordingly, finite am¬ 
plitude perturbations are required to trigger turbulence, 
a scenario referred to as by-pass transition 0. Many 
studies in the above flows have shown that the transition 
to turbulence is associated with the presence of 3-d exact 
coherent states Q. They appear in saddle-node bifurca¬ 
tions which in the state space of the system create regions 
of initial conditions that do not decay to the laminar pro¬ 
file, but instead are attracted towards the node-state . 
As the Reynolds number increases, the region widens, 
the node state undergoes further bifurcations and chaotic 
attractors or saddles are formed Initial conditions 

can only trigger turbulence when they reach into that in¬ 
terior region, i.e. cross the stable manifold of the saddle 
state on the boundary of the region Q. An “optimal” 
perturbation is one that can trigger turbulence and at 
the same time is a minimum of a prescribed functional. 
Popular is an optimization based on amplification or en¬ 
ergy gain over a given interval in time [9l4l7l| or on the 
total time-averaged dissipation 18, [l^. Because they 
take the time-evolution into account, they connect to op¬ 
timization problems in control theory BEl. 


We simplify matters here and focus on a geometric 
optimization by identifying initial conditions that will 
eventually become turbulent, without regard of the time 
it takes for them to become turbulent. The states are 
optimized so that a certain quadratic function, such as 
energy content or dissipation, is extremal: it is a maxi¬ 
mum in the sense that all initial conditions with a lower 
value of the quadratic function will not become turbu¬ 
lent, and it is minimal in that the first initial conditions 
that become turbulent have values larger than this opti¬ 
mum. At the optimal value there will then be at least one 
trajectory which neither becomes turbulent nor returns 
to the laminar profile: it lies on the stable manifold of 


the edge state [8|, so that the optimum is reached when 
the iso-contours of the optimization functional touch the 
stable manifold of the edge state (similar descri ptio ns of 
the state space structure can be found in [l3, 1^ 22. 2^1. 


THE MODEL 


To fix ideas and to keep the mathematics as simple 
as possible, we take the 2-d model introduced by Baggett 
and Trefethen 2^. The model we use is one of a set of 
many low-dimensional models of various levels of com¬ 
plexity 13, 22, iB-Slj. It has a non-normal linear part 
and an energy conserving nonlinearity, and, this being 
the most important feature for the present application, 
it is 2-d so that the entire phase space can be visual¬ 
ized (a property it shares with the illustrative model of 
[iiijll. Despite its simplicity, the model can be used to il¬ 
lustrate several features of the transition mechanisms in 
shear flows. 

The model has two variables, which may be thought 
of as measuring the amplitudes of streaks x and vortices 
y (see also [H), and one parameter R that plays the role 
of the Reynolds number: 


X = —x/R y — y^x^ -j- y^ 
y = -2y/R + x\/x^ P 


( 1 ) 

( 2 ) 


In order to highlight more clearly what happens near 
the origin, we magnify by rescaling the variables with 
the Reynolds number R (see [s^), i.e. we redefine the 
amplitudes x = x' jR?, y = y'/R? and the time t = Rt' 
such that (with the primes dropped) 

X = -X P Ry - y\Jx^ y^/R (3) 

y = -2y + x\Jx^ P y'^/R (4) 

Time evolution under the nonlinear terms alone preserves 
x^ Py'^, which may be thought of as a kind of energy, so 
that the nonlinear terms are “energy” conserving. For 
R < Rc = the only fixed point is x = y = 0, 
henceforth referred to as the “laminar” fixed point. At 
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Fig. 1. State space of the 2-d model for the transition to tur¬ 
bulence for R = 3. The open symbols mark the stable fixed 
point at the origin (“laminar” state) and the two nodes from 
the bifurcation (“turbulent” states). The full symbols are the 
edge states, and the red lines indicate the stable manifolds of 
the edge states. The black circle and the grey ellipsoid indi¬ 
cate the states where the energy 0 and the noise functional 
d are minimal, respectively. The points where they touch 
the stable manifolds are indicated by stars. 


R = Rc, symmetry related fixed points appear at {Xc, Uc) 
and {-Xc,-yc), with 


Xc = R{2R ± 2y/R^ - 8)/D± (5) 

= i?(i?2 _ 4 ± _ 8)(6) 


where D± = \/8 -I- 2R’^ ± 2i?V R^ — 8. The two fixed 
points closest to the origin are unstable, hence are sad¬ 
dle states, and the two further out are stable and hence 
nodes. The saddle states are the “edge states” Q and 
the node states are in the regions where turbulence would 
form, if more degrees of freedom where available. Never¬ 
theless, we will refer to them as the “turbulent” states. 

As i? —7^ oo, the saddles are to leading order in 1/i? 
located at ±(2,2/i?), which in the original coordinates 
represents an approach to the origin like ±(2/i?^, 2/i?^). 
The stable manifolds rotate so as to become parallel to 
the a;-axis, as we will see in the following. 


OPTIMAL INITIAL CONDITIONS OF MINIMAL 
ENERGY 

The Euclidean distance to the origin can be obtained 
from a quadratic form 

E={l/2)[x^ + y^) (7) 

which has the form of a kinetic energy. This assignment is 
further supported by the observation that E is preserved 
under time evolution by the nonlinear terms alone. In 
the sense described in the introduction, optimality with 
respect to this energy functional thus means the largest 
value up to which all trajectories return to the laminar 
state, and the smallest one where the first trajectories 
that evolve towards the turbulent state become possible. 
On the boundary between these two cases are states that 
neither return to laminar nor become turbulent, that lie 
on the stable manifold of the edge state. Geometrically, 
we are thus looking for the circle with largest radius that 
we can draw around the origin that just touches the sta¬ 
ble manifold. Algorithmically, we find this point by a 
modified edge tracking which minimizes the energy 0 
as described in the appendix. 

An example of such an optimal circle is given in Fig¬ 
ure 1, and its variation with R is shown in Figure 2. As 
the Reynolds number increases, the fixed point moves 
towards (2,0) on the abscissa, and the stable manifold 
rotates to being parallel to the abscissa. The point of 
contact between circles of equal energy and the stable 
manifold moves away from the edge state, approaches 
the y-axis and moves inwards to the origin like 1/i?. 

In an insightful discussion of the energy functional, 
Cossu [ 2 ^ notes that in the time-derivative of the energy 
functional only the linear parts of the equations of mo¬ 
tion remain and that the nonlinear ones drop out because 
energy is preserved. This observation allows to define 
a necessary condition for the location of the extremum, 
which for the two-dimensional example studied here im¬ 
plies that the optimum lies along the line connecting the 
laminar and the turbulent fixed points. One could then 
find the optimum by a one-dimensional search along this 
line. However, we did not pursue this further, as we also 
want to find optima with respect to other functionals that 
are not preserved by the nonlinear terms. 


OPTIMAL INITIAL CONDITIONS OF MINIMAL 
ENERGY DISSIPATION 

The diagonal terms in the linear part of the equations 
of motion correspond to the dissipation in the original 
Navier-Stokes equation. Accordingly, we can dehne a 
dissipation functional 0 

e=(l/2)(x2 + 2y2). 


( 8 ) 
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Fig. 2. Optimal states in energy for different R. The open 
symbol in the middle is the laminar state. States of fixed en¬ 
ergy are indicated by circles, and the points where they touch 
the stable manifolds (red lines) of the edge states (indicated 
as full symbols) are the points marked by stars. One notes 
that as R increases, the manifolds become more parallel to the 
i-axis, and the point of contact approaches the origin from 
the y-axis. 


Fig. 3. Optimal states with respect to the dissipation func¬ 
tional for different R. States of fixed dissipation are indicated 
by ellipses, and the points where they touch the stable mani¬ 
folds (red lines) of the edge states (full symbols) are marked 
by stars. The open symbol in the middle is the laminar state. 
One notes that as R increases, the point of contact moves 
very much as in the case of the energy functional. 


and study initial conditions that are minimal or optimal 
with respect to this functional. As in the previous exam¬ 
ple, the geometrical condition is that we now have to find 
the point where an ellipse touches the stable manifolds. 
This gives the ellipses shown for different R in Figure 3. 
Note that the points where the ellipses touch the stable 
manifolds are different from the ones of the energy func¬ 
tional, but their asymptotic behaviour for large R seems 
to be similar (see below). 


OPTIMAL NOISY TRANSITIONS 

As a third example we consider noise-driven transi¬ 
tions. To this end, the equations of motion are expanded 
to include a stochastic forcing of the individual terms, 

X = -X + Ry - yyj+ y'^/R + (9) 

y = -2y + x\/x'^ + y"^ / R + (10) 

where the noise is characterized by (^i) = 0 and 
= DijSij. We consider the case Du = 
D 22 = D, so that both components are driven with equal 
noise amplitude. In a linear approximation around the 
laminar fixed point, the non-normal coupling between the 
two components results in a probability density function 


(pdf) for the two components that is Gaussian with a 
covariance matrix given by 341 


p{x,y) = —[det((5)]2 exp{-QN{x,y)/D) 


( 11 ) 


with 


1/9 -3i? 

^ “ 2(i?2 + 9) \-SR 3A2 + 18 

and 

QN{x,y) = (13) 

= 2{R^+ 9 ) ~ ^ 6 ) j / 2 )( 14 ) 

Asymptotically, for i? —>■ 00 , the quadratic form becomes 
QN{x,y) —>■ 3y‘^f2, so that the Gaussian stretches out 
along the x-direction for increasing R. 

In the noisy case, transition is induced when a fluc¬ 
tuation carries the system across the stable manifold. 
A good estimate of the likelihood of transition can be 
obtained by considering the probability density at the 
transition point. Given the functional form of the pdf, 
the biggest contribution to its variations comes from the 
quadratic form in the exponent. The equation shows 
that the iso-contours p = const are ellipsoids determined 
by Qw = const that decrease or increase with the noise 
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Fig. 4. Probability density for the linearized equations with 
noise for iZ = 6 in the region of phase space where the tran¬ 
sition is expected to occur. For the statistics we calculated 
the time evolution of 20000 initial conditions starting at the 
laminar state for 20 time units with a step size of dt — 10“^. 
It can be seen that the iso-contours p = const are of elliptical 
shape. The star indicates the point where the noise functional 
touches the stable manifold of the edge state. 
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expected to occur. It is obtained by integrating 20000 
initial conditions in time for 20 time units starting at 
the laminar state with a step size of dt = 10“^. As the 
phase space is symmetric with respect to the origin, tra¬ 
jectories from the third quadrant are mirrored into the 
first quadrant. Figure |H obtained without the nonlin¬ 
ear part, shows the Gaussian shape of the iso-contours. 
Out of the 400 x 10® calculated points of the trajectories, 
more than 108 x 10® lay in the plotted region of phase 
space. In Figure [5] the nonlinear part is added and the 
pdf stretches out along the path to the turbulent state. 
The figure shows clearly that this happens close to the 
point where the ellipsoid QN{x,y) = const touches the 
stable manifold. Here more than 75 x 10® points lay in the 
interesting region of phase space. We note that the shape 
of the iso-contours of the pdf is independent of the noise 
amplitude D (within the linear approximation), so that 
changes in D will predominantly influence the likelihood 
of a transition, but not the path it takes. 

More examples of such iso-contours are shown in Fig¬ 
ures [T] and [S] for different values of R. With increasing 
R the ellipsoids become more elongated in x-direction, 
as a result of the asymptotic behaviour of the quadratic 
form noticed above (fT4l) . As they are rotated in the di¬ 
rection opposite to the rotation of the stable manifolds, 
the point of contact stays close to the edge state. Within 
the hydrodynamic interpretation, the transition is domi¬ 
nated by the streaks (x-component) that form as a result 
of the vortices (y-component), not by the vortices them¬ 
selves. 


SUMMARY AND CONCLUSIONS 


Fig. 5. Probability density as in Figure |4l but for the full 
nonlinear equations with noise. Note that the iso-contours 
of equal probability are stretched out towards the turbulent 
state and that they cross the stable manifold close to the point 
of contact indicated as the optimal state. 


amplitude D. Therefore, if we want to describe where 
a noisy trajectory crosses over to the turbulent state, 
we again have to study iso-contours of a quadratic form, 
Qn = const, and determine where they touch the stable 
manifold of the saddle state. In contrast to the energy 
functional © and the dissipation functional ([5|) , the fluc¬ 
tuation functional Qn depends on the Reynolds number. 
The point of contact between the probability iso-contours 
and the stable manifold then corresponds to the point 
where trajectories are most likely to cross over the stable 
manifold and to become turbulent. Alternatively, if one 
wants to push the system to become turbulent, small per¬ 
turbations in that region are most effective because the 
border is so close. 

Figures 0] and [5] show the relative probability den¬ 
sity to be at {x, y) in the region where the transition is 


The calculations illustrate how different optimization 
criteria select different optimal initial conditions for the 
transition to turbulence. Geometrically, this is to be ex¬ 
pected since different quadratic forms give rise to differ¬ 
ent ellipsoids in their iso-contours and hence also different 
points of contact with the stable manifolds. We note that 
the results of 13| suggest that for optimization with a 
time-integrated functional the difference between energy 
and dissipation functionals are smaller and may actually 
vanish. However, we have not pursued this question fur¬ 
ther. 

The variation of the optimal points of contact is sum¬ 
marized in Fig. [T] The data indicates that while the 
optimal perturbations are vortex like for the energy and 
the dissipation functional, they are streak like for the 
noisy transition. The difference can be rationalized by 
the different dynamics. In the deterministic cases, with 
the energy and the dissipation functional, small vortex 
like initial conditions can grow in time to develop the 
streaks which then drive the transition. The noisy sys¬ 
tem is always exposed to small perturbations which can 
grow to develop streaks, so that the pdf is elongated in 
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Fig. 6. Isocontours of the probability density function for dif¬ 
ferent R. Note that as R increases, the ellipsoids of the iso¬ 
contour become narrower in the y-direction and stretch out 
along the ^-direction. In combination with the rotation of the 
stable manifolds (red lines) the point of contact (stars) now 
stays close to the edge state (full symbols) and moves towards 
the x-axis. This is physically plausible, as a small perturba¬ 
tion in y (in the vortex direction) will produce a strong streak 
in the x-direction, and it is then the streak that triggers the 
transition. 


Fig. 7. Location of the points of contact. The top graph shows 
that for the energy and dissipation functional the optimal 
perturbations are vortex like and move towards the origin 
along the y-axis. The noise-optimals are more streak-like and 
move towards the x-axis, and actually remain close to the 
edge state. The edge state (red bullets) approaches (2, 0) as 
i? —>■ 00. 


the streak direction by non-normal amplification. There¬ 
fore, the transition happens on top of the already existing 
streaks and noise driven flows 351-1371 may show different 
structures at the point of transition than flows driven by 
judiciously chosen initial conditions. 

A final quantity to study is the scaling of the func¬ 
tional with Reynolds number, as shown in Fig. |51 De¬ 
spite the differences in dynamics, the functionals scale in 
all three cases like 1/i?^ for large R. The particular ex¬ 
ponent is specific to the model studied here and the type 
and form of the nonlinear interactions, as other nonlin¬ 
earities can require a rescaling near the origin 3^ . How¬ 


ever, the fact that all three cases show the same scaling 
could also apply to the full flow cases, as it is a conse¬ 
quence of the measure used and not the particular non¬ 
linearity at play. What the model also shows is that 
deviations from the asymptotic behavior appear close to 
the point of bifurcation. It is tempting to speculate that 
such effects may be responsible for the different critical 
exponents that have been observed in pipe flow or plane 
Couette flow, but that clearly requires the transfer of the 
present analysis to realistic flow simulations and a careful 
analysis of the asymptotic properties. 

The analysis of simple models has repeatedly helped 



Fig. 8. Scaling of the functionals for the optimal perturba¬ 
tions. Even though the point of contact moves differently, the 
critical values decay like \/R? in all three cases. 
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to elucidate many features of the transition to turbu¬ 
lence in shear flows, and to develop tools to explore them 
13) 22, 25- ^ I . It is in this spirit that we have used a for¬ 


ward integration technique to find the optimal points on 
the stable manifolds for different functionals and to ex¬ 
plore the changes with Reynolds number. We expect that 
many of the features described here can also be found in 
the high-dimensional state spaces of realistic shear flows, 
perhaps after suitable modifications and adaptations of 
the methods used to explore the high-dimensional spaces. 


APPENDIX 

In this appendix we discuss the modification of the 
edge tracking algorithm @ used for the determination of 
the initial conditions on the edge that optimize a pre¬ 
scribed quadratic functional Qpf{x,y). The functional 
may be the energy ©) the dissipation (|S]) or the argu¬ 
ment in the pdf (1121) . To keep the notation compact, we 
denote the equations of motion in vectorial notation as 
k = f. 

We begin with an arbitrary initial condition in the 
vicinity of the edge and we let it evolve in time to¬ 
wards the edge state. Unlike other edge tracking meth¬ 
ods, where trajectories are integrated until they are suf¬ 
ficiently close to the laminar or the turbulent state, we 
here stop the integration at the time when the distance 
to the edge state is minimal. The trajectory’s velocity at 
the turning point is then projected onto the normal of the 
stable eigenvector to decide if the tested initial condition 
moves upwards or downwards, towards the turbulent or 
the laminar state. With this criterion we can determine 
a point xq on the edge also when the point is very close 
to the edge and the time needed to pass the edge state 
becomes excessively large. 

We then propagate this point along the time direc¬ 
tion, xi = Xq -b sf(xo) by an amount s that is chosen 
such that Q{s) is minimized. Formally, 

_ F(xo)Qxo- bxg^Qf(xo) 

^ 2f^(xo) (5f(xo) 

In numerical implementations, ||sf|| is kept below a cer¬ 
tain threshold to stay in a region where linear approxi¬ 
mations are possible. Then a new edge tracking is started 
from xi and the process is repeated until the norm of the 
total shift ||sf|| falls below a convergence threshold, here 
10 -^ 
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